Criticality of a dissipative self-organized branching process in a dynamic population 
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We derive a general formulation of the self-organized branching process by considering sandpile 
dynamics in an evolving population characterized by "birth" (excitation) and "death" (de-excitation) 
of active sites (z = 1). New active sites are born in empty sites (z = 0) with a probability of n, 
whereas active sites die, thus becoming empty, with a probability A. Subsequently, when an active 
site becomes unstable (z = 2), it topples by transferring two grains to two randomly chosen sites 
with probability a or, by transferring only one grain to a randomly selected site (while retaining 
the other) with probability /3 = 1 + — — 2a, thus remaining active after toppling. We show that 
when sandpile dynamics occurs in an evolving population, self-organized criticality, characterized by 
a power-law avalanche size distribution with exponent t 3 — 3/2 and power-law avalanche duration 
distribution with exponent tt — 2 at very high dimension n » 1, is achieved even in the presence 
of dissipation (e = 1 — a — /3 > 0), contrary to previous claims. 
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The study of self-organized criticality (SOC) has re- 
ceived widespread attention because of its inherent sim- 
plicity of simulating scaling phenomena in a multitude of 
real systems. The sandpile model of Bak et al. is an im- 
plementation of SOC. It had a significant impact to the 
research community because of its capability of simulat- 
ing power-law distributions observed in a lot of systems 
in nature. How exactly the system approaches the SOC 
state is not so clear. A mean-field approach sheds light 
on this issue by providing a qualitative description of the 
spontaneous emergence of criticality in the system. Re- 
cently, Lauritsen, Zapperi and Stanley [l| proposed an 
analytic method that treats the sandpile dynamics as a 
branching process - a well-founded mathematical theory 
that was initially applied to understanding spread of epi- 
demics 2]. The model was aptly named self-organized 
branching process or SOBP, which is a very useful way of 
analyzing the sandpile dynamics. 

The SOBP relaxation process is conservative such that 
the transferred "energy" (grain) is equal to the "energy" 
lost by a relaxing site. Consequently, the system self- 
organizes to a critical steady state that is characterized 
by a power-law distribution for avalanche sizes (with ex- 
ponent r s = 3/2) and for avalanche lifetime distribution 
(with exponent tt = 2). The same group later demon- 
strated that any degree of dissipation in the relaxation 
rules leads to a subcritical steady state that manifests 
with a characteristic size in the avalanche distribution 
(power law followed by a significant drop-off after some 
characteristic avalanche size) [3j. A probability e that 
an unstable active site relaxes by dissipating grains out 
of the system was defined. Only at e = is the self- 
organized branching process become critical. The results 
apparently suggest that SOC can only occur in a conser- 
vative wherein the relaxation process is non-dissipative. 
Therefore, the existence of SOC in nonconservative mod- 
els would be highly desirable in this context, since, in 



practice, some type of dissipation is always present in 
nature gj|. An actual example of a dissipative system 
displaying SOC is an earthquake, which is successfully 
modelled by the Olami-Feder-Christensen (OFC) model. 
The OFC model sparked a wide interest because the re- 
laxation process is dissipative. However, the criticality 
as observed in the power-law earthquake intensity dis- 
tribution was recently questioned by de Carvalho and 
Prado |j|, who strongly claimed that the OFC model 
displays true criticality only for a non-dissipative relax- 
ation process. The apparent scale invariance observed in 
the avalanche size distributions produced from numerical 
simulations of the model is due to very long (but finite) 
correlation length. True criticality is supposed to display 
an infinitely long correlation length. Thus, the issue of 
observing SOC in dissipative systems remains an open 
question. 

In this Letter, we derive a general SOBP formulation 
that is carried out by first generalizing the toppling rules: 
An active site (z = 1) may transfer two grains to two 
randomly selected neighbors with probability a or, with 
probability [3, may transfer only one grain to a randomly 
chosen site while retaining the other, thus remaining ac- 
tive after toppling. Most importantly, we show that crit- 
icality is preserved, even with the presence of dissipation 
(a + (3 < 1), in a dynamically evolving population char- 
acterized by a birth rate 77 of new active sites and mortal- 
ity rate A of existing active sites, in contrast to previous 
claims |3j • The criticality of the sandpile is determined 
by the demographic parameters rather than the a nor (3. 

Demographically dynamic SOBP. — The mean-field de- 
scription of the sandpile model is defined by two stable 
heights, z — (inactive or empty state) and z — \ (ac- 
tive or occupied state). An active site becomes unstable 
by the addition of one grain, i.e., z = 2. In a similar 
fashion as in [lj, boundary conditions are implemented 



by allowing only n generations of toppling. This allows 
the mapping of the sandpile model to a branching tree 
with N = 2 n+1 - 1 nodes. 

Let us define p — Nj_/N as the density of active (z = 1) 
sites, where Aq is the actual number of active sites and 
N is the total number of sites. This density is time- 
dependent since Aq changes in time, starting at Aq = 
without loss of generality. Each of the (1 — p)N inactive 
(z = 0) sites has a probability A of being excited, whereas 
each of the pN active sites has a probability n of being de- 
excited. With such process, the system can be described 
as demographically dynamic. This update rule initiates 
each time step. 

Subsequently, an external energy grain is added to a 
randomly selected site. This grain has a probability p of 
landing on an active site [z = 1) thereby rendering it un- 
stable (z — 2). When this occurs, the demographic pro- 
cess is momentarily frozen and relaxation ensues. There 
are two ways that an unstable active site may relax: (1) 
with probability a the site topples by transferring two 
grains to two randomly selected neighbors, and conse- 
quently becomes inactive; (2) with probability /3 the site 
topples by transferring only one grain to a randomly se- 
lected neighbor, and remains active after toppling. The 
topplings result to an avalanche which propagates until 
all the sites relax to stable states, after which the demo- 
graphic process resumes. 

Since neighbor selection is uniformly random, the 
avalanche can be approximated as a branching process |2j 
with branching probability q^ = a5k,2 + P$k,x + (I — a — 
j3)8kfi- It follows that the generating function for the 
cluster size is 
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where a = (3 2 p 2 — 4ap (1 — (a + /3)p) and b = 0p. The 
branching process is critical 2] if (k) = ^2 kqk — 1 which 
corresponds to the condition p c = (2a + /3) _1 . It follows 
that this condition leads to a = 2b — 1. 

The avalanche size distribution P(s) is related to the 
generating function according to: 



By expanding Q around u> — and then comparing it 
with J2J, we are able to derive a recurrence relation for 
the cluster size distribution subject to the end conditions 
P(0) = and P(l) = (b 2 - a) /4a: 



= b(2s-l)P( s -l)-a( s -2)P( s -2) _ 



Graphically, Equation J2J) can be shown to define a 
power law with exponent r s = 3/2 for s >> 1 at the 
critical condition a — 2b— 1 with b G [0, 1). In particular, 
for b = the recurrence relation reduces to the analytic 
prediction in 0,|a], which only holds for odd values of s. 

Dissipation in the model is implemented by allowing 
Qk=o > or equivalently, a + [3 < 1. Hence, with prob- 
ability q = 1 — (a + /3), an unstable active site (z = 2) 
dissipates two grains out of the system. In [3j the dissi- 
pation probability qo = e. For conservative models, the 
amount of transferred grains is equal to the number of 
grains lost by a relaxing site, at least on the average jo], 
so that e = 0. 

Since we have embedded avalanche dynamics within 
a "demographic" process, then we define the relation 
X/rj = 2a + /3 — 1, which describes how avalanche dy- 
namics is interspersed in a dynamic population. The 
demographic parameters rj and A, and the branching pa- 
rameters a and (3 become interdependent. This condition 
can be equivalently written as, 
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by recognizing that p c — (2a + /9) _1 . Letting A — > and 
77 — > simultaneously is equivalent to the limit 1 — p c — > 
p c which implies that p c — > 1/2. This is exactly the result 
presented in [l| which can only be achieved if a — 1 and 
(3 = 0. For this case, the stationary state is a SOC state 
only in the non-dissipative regime (a + (3 = 1) as proven 
in 0. 

The rate of change of the density p of active sites can 
be written as 
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for which, cr(p, t)/N is a noise term corresponding to the 
fluctuation around the stationary value of p that van- 
ishes with increasing system size iV [3j. In order for the 



sandpile dynamics to be self-organized critical (SOC), 
the stationary state of the dynamical system must be a 
critical state. To show this, we solve for the fixed point 



solution of Equation (0 . Generally, the fixed-point solu- 
tion of the density of active sites p* does not coincide 
with its critical value p c . Moreover, p* < p c implies 
that the steady state is a subritical state for which the 
avalanche distribution possesses a characteristic scale and 
finite correlation length. Shown in Figgis a comparison 
of graph for dp/dt and the critical conditions a — 2b + 1 
and 1 — (2a + 0)p as functions of p. If the steady state 
is a critical state, all three functions must intersect at 
zero at p c . For a static population [3j, dp/dt = has 
a root p* < p c unless e = 1 — a — (3 — 0. This im- 
plies that self-organized criticality is only attained for a 
non-dissipative (e = 0) system as suggested in |3J. How- 
ever, in a dynamic population the critical state for which 
p c = (2a + /?) _1 is reached above a threshold excitation 
rate r\th that approaches to zero as the system dimen- 
sion n — ► oo. Illustrated in Fig. [3 is a transition from 
a subcritical steady-state to a critical steady-state for 
a degree of dissipation e = 0.2 for which a = 0.3 and 
fi = 0.5. With the coupling defined in Equation Q), 
the mean number of new sites gaining activity becomes 
equal to the mean number of active sites losing activity. 
Hence, the critical state coincides with a stationary state, 
that is, (1 — p)i] — pX — in Equation (JjjJ after carry- 
ing out a substitution using Equation fijl. Conversely, 
for a stationary state of the system wherein birth of ac- 
tivity equals death of activity, the SOBP can manifest 
criticality even in the presence of dissipation since a + /3 
can possibly be less than 1. The presence of a thresh- 
old excitation rate, rjth, is a requirement of time-scale 
separation between the avalanche propagation and the 
population dynamics ||j . We demonstrate our prediction 
with numerical simulations. The results are illustrated in 
Figure clearly showing the avalanche size and duration 
distributions with exponents t s — 3/2 and tt = 2 even 
with the presence of dissipation. 

Hence, a demographically dynamic population is capa- 
ble of achieving the critical state amidst the presence of 
local dissipation. This is significant for the fact that crit- 
icality has been observed ubiquitously in systems which 
do not display any degree of local conservation. More- 
over, these systems possess population dynamics. Several 
attempts have been proposed to show nonconservative 
systems exhibiting self-organized criticality. One intrigu- 
ing example is the Olami-Feder-Christensen model [?J, 
which quite recently had been proven to display critical- 
ity only in the conservative regime jj| . 

The idea of coupling population dynamics with self- 
organized criticality is reminiscent of the forest-fire 
model. But unlike the forest-fire model, the driving rates 
r] and A do not have to satisfy the limits A — » 0, r\ — ► 0, 
and X/rj — ► in order to display criticality in the sys- 
tem. As long as these demographic parameters satisfy 
a matching condition (i.e., X/rj = (1 — p c )/p c , wherein 
p c = (2a + /3) _1 , criticality is possible to achieve for 
V > nth- 



Fish school size scaling. — Bonabeau et al. |8[ presented 
data of aggregation of fishes and formulated an elemen- 
tary model of animal aggregation. In their mean-field 
model, they considered an interplay between aggregation 
and splitting. These parameters are analogous to the ex- 
citation rate ij and death rate A. In Figure 0J the data 
for tuna fish schools is fitted with the mean field solution 
(Equation J2J). The school size distribution manifests 
subcriticality owing to the mismatch between the demo- 
graphic parameters and the branching parameters. 



Adding a grain to a sandpile is analogous to stimu- 
lus reception of an ensemble of fishes. When a fish (or 
a small group of fishes) from the ensemble is excited by 
the stimulus it passes the information to other fishes. An 
information cascade ensues eventually synchronizing the 
action of the fishes via a proposed clustering biomech- 
anism known as allelomimesis 9] from which emerges 
aggregation to a food source as an example. Allelomime- 
sis has been found to be an appropriate mechanism that 
produces scaling in a wide sample of animal and human 
group size distributions |l(J . The subsequent stimulation 
of other fishes due to the observance of prior stimulated 
fishes can be thought of as a branching process wherein 
a stimulated fish has an average probability q^ of stim- 
ulating k other fishes (at random). With a probability 
<Zfc=o > 0, the relaxation is dissipative 3]. Unstimulated 
fishes (z = 0) of the ensemble may spontaneously get 
stimulated at an average rate of r\. On the other hand, 
stimulated fishes (z = 1) may spontaneously get deex- 
cited at an average rate of A. The dissipative nature of 
fish clustering is due to an imperfect tendency to copy 
neighbors as shown by Juanico et al. [3, ll3 through an 
agent-based model of clustering by allelomimesis. 



The resulting cluster size distribution was found to be- 
have as a power-law with an exponent 3/2. However, 
actual data exhibits a characteristic scale in the distribu- 
tion (subcriticality) that manifest as a cutoff size |8| ■ In 
Bonabeau et al.'s model, the truncation of the distribu- 
tion is attributed to a splitting parameter a which con- 
strains the lifetime of the fish school by virtue of breaking 
up after a certain amount of time. This splitting param- 
eter is analogous to deexcitation in our model that hap- 
pens with an average probability A. Subcriticality that 
manifest in the cluster size distribution is interpreted as 
the condition X/rj > (1 — p c )/p c . 



The main drawback with present SOC models is the 
need for some type of local conservation as an essential 
ingredient for a system to display SOC. The existence of 
SOC in nonconservative models would be highly desirable 
in this context, since, in practice, some type of dissipation 
is always present in nature |^. 
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FIG. 1: Rate of change dp/dt of the density of active sites 
(solid lines) , criticality condition 1 + a — 26 (dashed lines) and 
the function 1 — (2a + (3)p (dotted lines) as a function of p. 
(a) a — 0.3, /3 = 0.5, r) — 0, (static population); (b) a = 0.3, 
13 = 0.5, ri = 2" 4 , A = r\\ (c) a = 0.3, = 0.5, »j = 2" 4 , 
A = (2a + /3 - 1)7?; (d) a = 1.0, /3 = 0.0, r\ = (static 
population). The self-organized branching process is critical 
if all three curves mutually intersect at zero at p = p c , which 
occurs for (c) and (d). 
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FIG. 2: Steady state density of active sites p* as a function of 
excitation rate rj for different system dimension n. Parameter 
values: a = 0.55, /3 = 0. The steady state is a subcritical state 
below a threshold excitation rate rj t h(n), which decreases to 
zero with increasing system dimension n. Above r/th, p* — >• Pc- 
A phase transition from a subcritical state to a critical state 
occurs at rjth- 
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FIG. 3: Numerical simulations of the demographically dy- 
namic SOBP for system dimension n = 16 after an integration 
time of 2 17 iterations. Avalanche size distribution for a — 0.3, 
13 = 0.5 (dissipative) and a = 1, j3 = (non-dissipative) . 
Self-organized criticality is achieved for a dissipative sandpilc 
(a + /3 < 1). At criticality, the distribution follows a power 
law with exponent r s = 3/2 that is consistent with mean-field 
analysis. The inset graph shows the avalanche duration dis- 
tribution that fits to a power law with exponent tt — 2 at 
criticality. 
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FIG. 4: Normalized fish school size distribution data from 8] 
fitted with Eq. 0: a = 0.25, /3 = 0.5, \/r) = 0.25. The 
truncation at the upper tail of the distribution is evidence of 
subcriticality owing to significant probability of deexcitation 
(A = 0.25??). 



